Journal of Molecular Evolution

Language genes can be tentatively considered as a subset of cognitive genes, although they are often discussed separately. During the evolution of SNVs (single nucleotide variations) in cognition-related genes, do language genes and cognitive genes exhibit significantly different intensities of change at several key evolutionary moments--namely, the inflection points or derivative peak positions of similarity curves drawn from multi-SNV locus bases across samples? In this study, nine distance/similarity metrics (Bray-Curtis, Cosine, Pearson, Spearman, Hamming, Jaccard, Matching, Kulczynski, and Gower) were employed to analyze 413 samples from 11 taxonomic groups, targeting SNV loci in language/cognition-related genes (13,415 effective loci, approximately 400 loci per gene), with pp6 (Homo_sapiens.GRCh38) as the reference. For each method, sample similarities (defined as 1/(1+distance)) were independently sorted in ascending order to generate raw similarity scatterplots. Due to the large sample size and representativeness, the scatter density on the similarity curves was high, and no smoothing was applied. Derivative values were calculated from adjacent similarity differences to identify peaks of evolutionary rate change (top 10 peaks per method). Combined with functional annotations of 33 language/cognition-related genes, we quantified the difference scores and occurrence frequencies of the two gene categories at the peak positions. The results indicate that cognitive-related genes exhibit slightly higher occurrence frequencies in peak windows and higher average difference scores per gene than language genes. Comparative analysis of SNVs at the peak samples and their left-side windows revealed that at positions 381-382, all nine methods shared three intersecting mutation loci, involving language genes (NFXL1, SRGAP2, SRGAP2C); at positions 355-356, there was one intersecting mutation locus, involving a language gene (SRGAP2). This suggests that certain mutations in language genes may have played a distinctive role at critical junctures in the evolution of cognitive abilities.

Since the discovery of the genetic code, a frequently discussed question is whether the numbers of synonymous codons for the various amino acids are randomly distributed or were shaped by evolutionary constraints. In this study, we analyze for the standard as well as alternative genetic codes, the correlations between codon degeneracy and amino acid frequencies in proteins (neglecting differences in gene expression). To taking into account the effect of GC content, expected codon frequencies rather than codon multiplicity need to be considered. A strong correlation of these frequencies with amino acid abundance is revealed. Furthermore, we identify consistent patterns of over- and underrepresentation of amino acids across domains of cellular life as well as viruses. For example, the codons for glutamate, aspartate, lysine, and methionine are consistently overrepresented across domains, while cysteine, arginine, histidine, and proline are underrepresented. Subsequently, we discuss the role of biosynthesis costs of amino acids and other factors such as the order of amino acid recruitment in early evolution, exposure to oxidative stress, and sulphur availability. We hypothesize that in a first phase of evolution, the genetic code evolved in a way so as to comply with the different demands for amino acids. In a second phase, after the code was frozen, changes in demand led to deviations from the strong correlation between codon multiplicity and amino acid frequency by natural selection. This comprehensive analysis offers new insights into the interplay between genetic code structure and amino acid usage across all domains of life.

1. Plants produce a great chemodiversity, which is the diversity of specialized metabolites (SMs). These SMs are produced in complex metabolic pathways and play an important role in inter-species interactions. There are numerous hypotheses about the evolutionary processes which brought about and maintain chemodiversity. Some have been partially tested in lab and field studies. However, some of their assumptions and predictions are better tested by quantitative modeling, and so far no quantitative model has investigated the role of metabolic pathways. 2. To close this gap, we developed an individual-based model for metabolic pathway evolution. It models enzymes creating metabolites with various modifications. Enzymes undergo inheritance and mutation. We used the model to compare the screening and interaction diversity hypotheses. 3. The screening hypothesis predicts promiscuous enzymes, genetic drift, the presence of many non-beneficial metabolites, and high metabolite richness. The interaction diversity hypothesis predicts specialized enzymes, selection, the almost exclusive presence of beneficial metabolites, and situation- dependent metabolite richness. We found that the patterns predicted by the screening hypothesis did not occur, while those predicted by the interaction diversity hypothesis did. 4. This provides reason to favor the interaction diversity hypothesis over the screening hypothesis when connecting empirical results to their evolutionary context

This study aims to quantify the genetic similarity of different species (from fish to humans) to the human reference genome (pp6, Homo sapiens.GRCh38) based on the allele presence/absence patterns of 33 language/cognition related gene SNV loci, identify key breakpoints during evolution, and evaluate the enrichment of language and cognition genes at these breakpoints. We designed a similarity calculation method relying on binary features (four columns for A/T/C/G), adopted five difference/distance measures (Sorensen, Rogers, Nei, Reynolds, and Hellinger), and converted them into similarity values (1/(1+distance)). For each method, samples were independently ranked, the first derivative of similarity was computed, and the top 12 peaks were selected as candidate breakpoints. Results show that the similarity curves from the five methods are highly consistent (correlation coefficients >0.9), with major peaks concentrated at positions 355, 363, 381, 382, 390, 400, etc., where the corresponding samples are predominantly ancient hominins and primates. Furthermore, we defined 13 peak groups (starting positions 355-401). For each peak within a group, pairwise SNV differences between the peak apex sample and its immediate left neighbor were compared, and the intersection F_INTERSECTION (shared differential loci) was obtained. For each F_INTERSECTION, we calculated the proportions of language genes and cognition genes. In addition, we computed the differential sets between adjacent groups' F_INTERSECTION to trace the gradual emergence of new loci. In F_INTERSECTION, language genes accounted for an average of 59.5%, and cognition genes for an average of 62.9%. The proportion of language genes reached a peak at position 383 (61.2%), while cognition genes peaked at position 386 (64.9%). High frequency peak samples include c25, c27, and ja2, suggesting that language cognition genes may have undergone independent intensification during Eurasian evolution. Differential analysis between adjacent F_INTERSECTION revealed a stepwise acquisition of new loci from position 355 to 401, with three bursts of newly added loci along the entire evolutionary axis. This study provides a quantitative framework based on similarity curves, offers a novel molecular perspective for understanding the evolution of language and cognitive abilities, and highlights the potential importance of East Asian archaic hominins in the evolution of language cognition genes.

We present a genome assembly from a female Cistugo seabrae (Seabras wing-gland bat; Chiroptera; Cistugidae). The genome sequence is 1.9 gigabases (Gb) in span. The majority of the assembly is scaffolded into 25 chromosomal pseudomolecules, with the XX chromosomes assembled. The assembly has a contig N50 of 51.9 Mb and a scaffold N50 of 91.4 Mb. Species taxonomyEukaryota; Metazoa; Chordata; Craniata; Vertebrata; Euteleostomi; Mammalia; Eutheria; Laurasiatheria; Chiroptera; Yangochiroptera; Vespertilionoidea; Cistugidae; Cistugo; Cistugo seabrae Thomas, 1912 (Teeling et al., 2005; Bickham et al., 2004; Lack et al., 2010).

Model selection is a fundamental step in phylogenetic analysis that determines the best-fit model of sequence evolution for a given multiple sequence alignment. Popular model selection methods, such as ModelFinder, rely on statistical information criteria, such as the Bayesian Information Criterion (BIC) or the Akaike Information Criterion (AIC). However, these approaches are computationally expensive and the use of information criteria has been the subject of ongoing discussion. Recently, machine learning has emerged as a promising approach for phylogenetic model selection in both nucleotide and protein sequence analyses. ModelDetector is currently the only machine learning-based method for amino acid substitution model selection. However, because ModelDetector was trained on simulated data, it does not perform well on real datasets. Another limitation is that it does not support different rate heterogeneity across sites (RHAS) models. To overcome these limitations, we introduce ProtFinder, an efficient machine learning framework for protein model selection that predicts amino acid substitution models, RHAS models, and amino acid frequency models. To enable ProtFinder to work with real datasets, we employed a transfer learning strategy consisting of three stages: (1) initial training on large-scale simulated data, (2) joint training on both simulated and real data, and (3) final fine-tuning using real data only. Experimental results show that ProtFinder outperformed ModelDetector in amino acid substitution model selection. ProtFinder achieved comparable accuracy to the maximum likelihood method ModelFinder for substitution model selection on medium and large MSAs. It performs slightly better than ModelFinder in RHAS model selection and substantially outperforms it in amino acid frequency model determination. Notably, ProtFinder is up to 1,400 times faster than ModelFinder in terms of inference time, making it particularly suitable for medium and large datasets.

Gene duplication is an important evolutionary mechanism contributing to genomic innovation and adaptation. While the long term fate of duplicated genes has been extensively studied, immediate effects such as those resulting from increase in gene dosage remain understudied. Of particular interest is the question of how expression level relates to the optimum after a duplication, as it could fall short, exceed or directly reach the optimal expression. In this study, we investigated the relationship between expression and fitness for 11 genes previously reported to increase fitness when duplicated in Saccharomyces cerevisiae to compare the expression level achieved by duplication to the one maximizing fitness. While most of the studied genes show little fitness variation across the expression levels tested, we found that ESP1, which encodes a separase essential to anaphase initiation, exhibits a clear fitness optimum at an expression level above the native single-copy state. However, the expression level associated with gene duplication exceeds the optimum. We find that ESP1 overexpression is associated with a reduced cell size and a lower proportion of the cell cycle spent in G2 phase, consistent with its role in cell cycle progression. We propose that the duplication of ESP1 accelerates the cell cycle, resulting in smaller, faster-dividing cells. When measured in different environments, the fitness functions also exhibit similar shapes, with slight variation in the optimal expression. Our findings highlight how gene duplication can confer an immediate selective advantage through dosage increase while being non-optimal, reversing the directionality of natural selection on expression level after duplication.

Gene model for the ortholog of Density regulated protein (DENR) in the Apr. 2013 (BCM-HGSC Dpse_3.0/DpseGB3) Genome Assembly (GenBank Accession: GCA_000001765.2) of Drosophila pseudoobscura. This ortholog was characterized as part of a developing dataset to study the evolution of the Insulin/insulin-like growth factor signaling pathway (IIS) across the genus Drosophila using the Genomics Education Partnership gene annotation protocol for Course-based Undergraduate Research Experiences.

Symbiosis has lasting effects on the genomic evolution of interacting organisms. Shifting niches and selection intensities are well-documented drivers of adaptive evolution, but nonadaptive evolution is an equally important, if less conspicuous, dimension of symbiont genome evolution. Through a reduction in effective population size (Ne), symbiosis can reduce the efficacy of natural selection in symbionts; however, evidence of this "genomic syndrome" is derived almost entirely from bacteria. Whether it extends to eukaryotic symbioses with complex demography remains unclear. Here, lichen-forming algae are examined for signatures of a genome-wide reduction in the efficacy of selection arising from demographic shifts associated with symbiosis. To test this, four lichen-forming algal taxa were compared to their closest free-living relatives using complementary measures of (i) the strength of molecular evolution (dN/dS, K), (ii) codon use bias (ENC), and (iii) base composition at synonymous sites (GC3). Lichen-forming algae showed heterogeneous responses across all signatures of molecular evolution examined. The effect of lifestyle (lichen-forming vs. free-living) on dN/dS ratios and codon use bias was lineage specific. GC3 was uniformly reduced in lichen-forming taxa, but the underlying causes differed, and no genes showed uniformly intensified or relaxed selection across all four lichen-forming taxa. These results suggest that, while symbiosis often reshapes symbiont population genetics in ways that elevate the role of drift, lichen-forming algae do not show a uniform reduction in the efficacy of selection. Rather, the demographic effects and functional demands of each specific lichen symbiosis likely shape genome evolution in lineage-specific ways.

Across metazoa, morphogen gradients differentially regulate gene expression and activate a spatially distinct program to specify body axis development. The Drosophila gap gene network, initiated by maternal morphogen Bicoid, is one of the most well-studied systems. Several regulatory interactions act synergistically to produce distinct gap gene expression patterns along the anterior-posterior (AP) axis of the blastoderm stage Drosophila embryo to ensure the proper segmentation of the larval and, eventually, adult stage fly. Several mathematical models have been proposed to summarize the interconnectivity of gap gene regulatory elements and predict expression in mutant systems. However, these models have not successfully predicted the gap gene expression profile over the entire AP axis. Here, we present an Enhancer Activity-informed Gap GEne Regulatory network (EAGER) model that incorporates enhancer activity-driven differential regulation along the AP axis, which successfully summarizes gap gene expression patterns over the entire AP axis. We validated the predictions of EAGER on Kr mutants and performed a comprehensive parametric sensitivity and identifiability analysis to evaluate the robustness of the EAGER model fits and predictions. We also propose a reduced version of the model, rEAGER, which identifies a minimal set of regulatory interactions and successfully summarizes the gap gene interactions over the entire AP axis. Our results suggest that the expression driven by individual enhancers must be accounted for in models of developmental pattern formation.

Single-cell sequencing provides a rich source of data that can be used to model and estimate the evolution of somatic cells within individuals. Many common models for modeling DNA evolution in phylogenetics are based on continuous-time Markov processes that operate on aligned sites. For haploid DNA, such as found in mitochondria, there are many well established models for these processes but since nuclear DNA is diploid, containing copies of both the maternal and paternal genome, it is natural to develop continuous-time Markov processes adapted to diploid DNA. One possible model is that both the maternal and paternal evolve independently according to haploid models so that, in principle, DNA comes with a phase. i.e. a label that tells if the site was derived from the maternal or paternal sequence. However, this information may not be available in practice and only an unphased diploid genotype is available. That is, an unphased model groups both of the possible maternal-paternal labels into a single state. This grouping process, called ``lumping'', is used by other methods as well and results in a coarser process. Care must be taken so that the coarser process can still be modeled as a continuous-time Markov chain; the lumped process of the original phased Markov chain may or may not be Markovian, depending on restrictions imposed on the generator matrix. Here, we explore the mathematical formulations of phased and unphased Markov processes, as well as a collection of other diploid DNA models from the literature.

Cells and tissues have the ability to respond to different stresses and maintain their homeostasis. The b-Zip transcription factor Xrp1 was found to play a pivotal role in this context in Drosophila. Xrp1 is required for the response to irradiation, in genetic conditions of ribosomal protein (Rp) deficiency leading to cell competition, and more generally in mediating the Integrated Stress Response (ISR). However, how these diverse stresses converge to and activate Xrp1 is still unclear. Our results demonstrate the existence of three separate mechanisms of Xrp1 induction. In response to Rp stress, a short isoform of Xrp1 is produced and requires the function of the ribosomal protein RpS12. In response to irradiation, a specific mRNA isoform is transcriptionally induced, which encodes the Xrp1-short protein isoform, in a p53-dependent manner. In response to the ISR, the long isoform of Xrp1 is translated by the derepression of a uORF present on the Xrp1 5'UTR, a mechanism that requires eIF2 phosphorylation, RpS12, and a transcriptional induction by ATF4. Therefore, our work provides a comprehensive understanding of the role of Xrp1 in different situations of cellular stress, illustrating the complex level of integration of the mechanisms controlling this central gene.

THAP9 is a transposable element-derived gene which encodes a protein that is homologous to the active Drosophila P-element transposase (DmTNP). Both THAP9 and DmTNP possess a C-terminal domain (CTD) which is functionally uncharacterized. Sequence and structural analysis suggest that the THAP9-CTD has a novel fold which is only found in THAP9 homologs. To explore the evolutionary history and characteristics of this novel domain, exhaustive phylogenetic analysis (using MSA, structure prediction, MSTA-based clustering) was performed. THAP9-CTD homologs were more widely distributed throughout the animal kingdom in comparison to DmTNP-CTD homologs which were restricted to arthropods. Moreover, the THAP9-CTD homologs were more conserved, especially among mammals and birds and their average length increased in a class-specific manner. Comparison with the DmTNP-CTD homologs demonstrates that although their respective CTDs may have evolved independently, they both surprisingly share similar secondary structure elements consisting of three conserved helical regions made of hydrophobic residues that are predicted to make up a conserved core. The role of the respective CTDs were further investigated by creating truncation mutants lacking the CTD. Interestingly both THAP9 and DmTNP truncation mutants are still capable of DNA excision and integration suggesting that their respective CTDs are not essential for DNA transposition. Moreover, CTD truncation favours DNA integration in THAP9: this suggests that CTD acquisition during evolution may have led to THAP9 domestication as observed in other transposable element-derived genes like Rag1 and piggybac, which have similar terminal regulatory domains.

The majority of eukaryotes encode some intron-containing pre-tRNAs. Splicing of these pre-tRNAs requires a dedicated tRNA splicing machinery. The fungal and trypanosome tRNA ligase, Trl1, and the human RNA ligase, RTCB, catalyze an essential step in tRNA splicing. However, Trl1 and RTCB are nonhomologous and biochemically and structurally distinct from each other. Therefore, Trl1 could serve as a broad-spectrum antifungal and anti-trypanosomal target. While the functions and requirements of the three catalytic Trl1 domains have been extensively characterized in the model yeast Saccharomyces cerevisiae, the roles of Trl1 orthologs in pathogenic fungi remain unexplored. Here, we validate Trl1 as one of the few promising novel drug targets for the development of antifungal therapeutics. Functional analyses of the three Trl1 domains show that only the "sealing" domain is essential for growth and viability in Candida albicans and Aspergillus fumigatus. In contrast, the two "healing" domains are dispensable in these pathogenic fungi, suggesting the presence of redundant healing enzymes, unlike in S. cerevisiae. These findings indicate that only the sealing domain is a good drug target. Our analysis also shows that the Mucor enzyme, which only contains the sealing domain, is essential. Using a Caenorhabditis elegans infection model of C. albicans, we further demonstrated that inhibiting Trl1 expression protects worms during an established infection. In contrast to these fungal pathogens, we show that all three domains of Trl1 are essential in Trypanosoma brucei. Our findings show that the essentiality of the Trl1 sealing is conserved in important human pathogens and provides an impetus for future drug development. SIGNIFICANCEFungal infections are an important cause of human disease and death and difficult to treat and there is an urgent need to develop additional drugs. Based on studies in yeast, one promising target for antifungal drug development is the tRNA splicing pathway. Human tRNA ligase is fundamentally distinct from the fungal one. To investigate the possibility of developing tRNA ligase-targeting drugs, we investigated the function of the catalytic domains of fungal tRNA ligase in different fungal pathogens. Surprisingly, only the first domain is essential in these pathogens and yeast is not a good model fungus. In contrast, all three domains of Trypanosome tRNA ligase are essential. These findings provide an impetus for future drug development.

Gene model for the ortholog of Insulin-like peptide 3 (Ilp3) in the D. pseudoobscura Apr. 2013 (BCM-HGSC Dpse_3.0/DpseGB3) Genome Assembly (GenBank Accession: GCA_000001765.2) of Drosophila pseudoobscura. This ortholog was characterized as part of a developing dataset to study the evolution of the Insulin/insulin-like growth factor signaling pathway (IIS) across the genus Drosophila using the Genomics Education Partnership gene annotation protocol for Course-based Undergraduate Research Experiences.

Two rounds of whole-genome duplication (2R-WGD) in the vertebrate lineage are thought to have driven gene diversification. Mammals have four genes encoding five nebulin superfamily proteins (LASP1, LASP2/NEBL, NEB, and NRAP) hypothesised to have arisen from invertebrate LASP via WGD. Here, we examined the evolutionary fate of the nebulin superfamily after an additional WGD, teleost-specific genome duplication (TGD). Genome and transcriptome analysis of zebrafish (Danio rerio) and medaka (Oryzias latipes), of which the complete genomes were recently determined, showed only lasp family genes retained duplicate copies. Zebrafish retained two LASP2/NEBL-derived ohnologs; one (LASP2L) had been misidentified as a nebulette gene in the Zebrafish Information Network (ZFIN) despite lacking nebulette-specific exons. Another (LASP2/NEBL) is newly identified in this study because it lies in a newly determined region of the latest genome assembly. Contrastingly, medaka retained two LASP1 ohnologs and one LASP2/NEBL. The duplicated lasp family proteins showed overlapping and specific tissue expression patterns. We confirmed all three zebrafish lasp proteins bound actin and co-localized with actin-rich structures when ectopically expressed in cultured fibroblasts. This supports that TGD-derived lasp duplicates underwent lineage-specific retention and expression divergence in teleosts whilst conserving core actin-binding properties.

In the genus Paramecium, the macronuclear genome is remarkably compact and optimized for gene expression. As a means to explore eukaryotic transcription in the context of a streamlined genome and shed light on the role of sequence architecture on gene expression and loss, we analyzed the distribution and diversity of candidate transcription initiation sites (TISs) in Paramecium sexaurelia, Paramecium tetraurelia and their outgroup, Paramecium caudatum. Our analysis suggests that for Paramecium, most genes have very short 5 prime UTRs (40 bp or less) and their transcription initiation regions (TIRs) have a median dispersion (akin to width) of 8-10 bp. The TIRs for the three species have high AT content. TIR dispersion is not to gene expression. However, mean TIS position relative to the translation start site per gene does in gene expression for the three species, and is often conserved between them. While mean TIS position and gene expression are linked, gene expression itself is the main driver of paralog retention in the aurelias. As compared to other eukaryotes, Paramecium has a uniquely well-defined and short main TIS region, and sequence motifs that likely diverge from the consensus in multicellular eukaryotes.

Forms I&prime, I&Prime, and I RuBisCOs were recently identified in bacterial metagenome-assembled genomes (MAGs) and are significant for evolutionarily bridging form I, which is responsible for the major part of carbon fixation on Earth, and other forms. Pioneering studies focused mainly on the biochemical and structural characteristics of "form I-related" RuBisCOs and deduced the evolution of this protein family over phylogenetic relationships among form I and form I-related sequences. We here surveyed form I-related sequences in public databases and conducted phylogenetic analyses incorporating the updated diversity of form I-related RuBisCOs. This study yielded two major outcomes regarding RuBisCO evolution. First, the analyses presented here suggest uncertainties overlooked in the previously published studies. It appeared to be difficult to conclude (i) the phylogenetic relationship among forms I, I&prime, and I&Prime, (ii) the monophyly of form I&prime, or (iii) the monophyly of form I. Second, we identified a gene encoding form I RuBisCO in an archaeal MAG, and this protein is likely involved in the pentose bisphosphate pathway. The results presented here suggest that the evolutionary route leading to form I in RuBisCO evolution remains unsettled.

Understanding the genomic consequences of thermal adaptation in fungi is crucial, as rising global temperatures are expected to have negative impacts on food safety and human health. The family Chaetomiaceae contains a large number of thermophilic fungal taxa, but previous studies have reported inconsistent optimal growth temperatures (OGT) for the same strains, obtained with various laboratory methods. Here we applied a standardized laboratory approach to measure OGT across strains of 17 Chaetomiaceae species and used a phylogenomic approach to test associations between OGT, rates of genome evolution, and strength of purifying selection. Compared to mesophiles, thermophilic fungi showed faster nucleotide substitution rates. In addition, thermophiles showed lower dN/dS ratios than mesophiles, suggesting stronger purifying selection on conserved orthologs. We hypothesize that the elevated substitution rates are linked to high growth rates, as thermophilic fungi grew significantly faster than mesophilic ones. Our results show that selective pressures may act at different temperatures for distinct genomic characteristics. Genome size was lower at OGT [≥] 35{degrees}C compared to mesophilic species, while GC content did not show a large difference between mesophiles and thermotolerant species, but increased in thermophilic species with an OGT [≥] 45{degrees}C.

The production of gametes is one of the universal processes of life and given millennia of evolution many of the core genes involved are expected to be highly refined and resistant to change. Here we examine the pattern of expression and functionally characterize three genes widely involved in gametogenesis: Vasa, Piwi, and Dnmt1. We examine their cellular location during both oogenesis and spermatogenesis in Oncopeltus fasciatus, a hemipteran insect, using fluorescent in situ hybridization chain reaction. In females, expression of Piwi was restricted to the trophocytes, but Vasa and Dnmt1 were expressed in developing ooctyes as well as the trophoctyes. In males, Piwi was expressed in the testis germline stem cells (GSCs), but Vasa and Dnmt1 expression was absent from these cells. Additionally, Piwi and Vasa were expressed in secondary spermatogonia and spermatocytes while Dnmt1 was predominantly expressed in the primary spermatogonia. We also functionally characterized their necessity for gametogenesis after RNAi-mediated gene expression knockdown. Vasa was not required for oogenesis but was required during spermatogenesis. Piwi and Dnmt1 were required for both oogenesis and spermatogenesis. These results were in contrast to their requirement for these processes from other organisms, which highlight the unexpected variation found in the processes and suggest that conservation depends on the level at which genes are examined: sequence, expression, or phenotypic and biochemical function. This suggests there is a more nuanced evolutionary story of conserved, yet plastic, gametogenic gene set in Metazoa that deserves further study across more organisms.